function elements = InitializeMaterial(elements,model)
ne = model.countElements();
for i=1:ne
    e = model.getElement(i-1);
    ng = e.getNumGaussPoint();
    for j=1:ng
        elements(i).constitutiveMatrix(:,:,j) = ...
            C_plane_strain(elements(i).E,elements(i).nu);
        elements(i).alpha(j) = 0;
        elements(i).alpha_old(j) = 0;
        elements(i).beta(:,j) = [0 0 0 0]';
        elements(i).beta_old(:,j) = [0 0 0 0]';
        elements(i).ep(:,j) = [0 0 0 0]';
        elements(i).ep_old(:,j) = [0 0 0 0]';
        elements(i).stress(:,j) = [0 0 0]';
    end
end
end
